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ABSTRACT 

X-ray flares have routinely been observed from the supermassive black hole, Sagittarius A* (Sgr A*), 
at our Galactic center. The nature of these flares remains largely unclear, despite of many theoretical 
models. In this paper, we study the statistical properties of the Sgr A* X-ray flares, by fitting the count 
rate (CR) distribution and the structure function (SF) of the light curve with a Markov Chain Monte 
Carlo (MCMC) method. With the 3 million second Chandra observations accumulated in the Sgr A* 

X-ray Visionary Project, we construct the theoretical light curves through Monte Carlo simulations. 

We find that the 2 — 8 keV X-ray light curve can be decomposed into a quiescent component with 
a constant count rate of ~ 6 x 1CD 3 count s _1 and a flare component with a power-law fluence 
distribution dN/dE oc E~ aE with cue — 1-65 ± 0.17. The duration-fluence correlation can also be 
modelled as a power-law T oc £“ ET with get < 0.55 (95% confidence). These statistical properties 
are consistent with the theoretical prediction of the self-organized criticality (SOC) system with the 
spatial dimension S = 3. We suggest that the X-ray flares represent plasmoid ejections driven by 
magnetic reconnection (similar to solar flares) in the accretion flow onto the black hole. 

Subject headings: Galaxy: center — methods: statistical black hole physics — accretion — accretion 
disks 


1. INTRODUCTION 

Sagittarius A* (Sgr A*) at the center of the Milky Way 
is an excellent laboratory for studying the accretion and 
ejection of matter by supermassive black holes (SMBHs). 
There have been quite a number of observational and the¬ 
oretic al s tudies of Sgr A* (see r eviews by Genzel et al. 
2010 and Yuan & Narayan 2014). The bolometric iumi- 


nosity of Sgr A* is Lboi ~ 10~ :j I/Edd (where LecM is the 
Eddington luminosity), which is five orders of magnitude 
lower than that predicted by a st andard thin disk accre - 
tion at the Bondi accretion rate (Baganoff et al.| 

We now understand that an advection-dominated 
tion flow scenario works for Sgr A*, and that 


2003). 


accre- 
tlie low 

luminosity is due to the combination of t he low radiat ive 
efficiency and the mass l oss via outflow ([Yuan fe Narayan 
2014| [Yuan et al.||2003[ |2012 |Narayan et al. 2012 ' 

al.||2013~|Wang et al.|2013' 


Li et 

__ Sgr A* is usually in a quies¬ 
cent state, and occasionally shows rapid flares (on time 
scales ~ 1 hour), most significantly in X-ray (Baganoff 


et al. 2001 1 and near-infrared (NIR; |Genzel et aT 2UU5f 


Chez et ai. 12004). The flare rate is roughly once per day 


in X-ray and more frequently in NIR. 

Theoretical interpretations of the flares include the 
electron acceleration by, e.g., shocks, magnetic reconnec- 


(e.g., 

Yuan et al. 2004 Eckart et al.|2006; Yusef-Zadeh et 

al. 2009 

Dodds-Eden et ai. 2009, 

2010; 

Chan et ai. 2015) 

or in an assumed jet (e.g., 

Markofi 

let ai. 

2001). The other 


scenarios in terms of the flaring cause are transient fea- 
tures in the accretion fl ow, such as accreti on instability 
dTagger fe Melia||2006|), orbiting hot spot (|Broderick fc 


expanding plasma blob (Yuset-Zadeh et alT 


20001 lEckart et al.|2006| Yusef-Zadeh et al.|2009||Trap et 


[20llf7and tidal disruption of asteroids (Zubovas et 

However, the nature of the flares is still under 


debate. 


detected (e.g., 

Baganoff et al. 2001 Porquet et al.| 2003 

Belanger et al. 

2005 

Eckart et al. 2006; Aharonian et ai. 

2008; Marrone et a 

. 2008; Porquet et al. 2008; Yusef- 

Zadeh et al. 2008 

jT\Z - 7TT — _ 1 V _ 

Legenaar et al. 2013j). Thanks to 


dfj ii 
er al. 


after XVB 1 in 2012, a total of 39 flares have been added 


(Neilsen et al. 2013). This substantially increased sam¬ 
ple size of the fl ares has enabl e d the statistical study of 


the properties. Neilsen et al. (2013) first analyzed the 


distributions of the fluences, durations, peak count rates 
and luminosities of the XVP detected flares. iNeilsen etl 
al. ([2015) further investigated the flux distribution of the 
XVP X-ray light curves taking into account the Possion 
fluctuations of the quiescent emission and the flare flux 
statistics. 

In this work we advance the statistical analysis of the 
XVP light curve by simultaneously using its flux distri- 


1 http://www.sgra-star.com/ 
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bution and structure function (SF). We generate X-ray 
light curves through Monte Carlo realizations of both 
the quiescent and flare contributions. The realization of 
flares assumes power law models for their fluence function 
and width-fluence dependence. The parameters of these 
power laws are in turn constrained by comparing the sim¬ 
ulated data with the observed one, using the Markov 
Chain Monte Carlo (MCMC) approach. 

The major motivation of the statistical analysis is 
to determine whether the SOC model can explain the 
Sgr A* X-ray flares, and if so to probe the dimensionality 
of the process leading to the flares. We compare our con¬ 
strained power law indexes (i.e., the power law indexes 
for the fluence distribution and the duration-fluence cor¬ 
relation) with the expectation of the self-organized criti¬ 
cality (SOC) theory, which describes a class of dynamical 
system with nonlinear energy dissipation that is slowly 
and continuously dri ven toward a critical value of an 


instability threshold (Katz 1986 Bak et al. 1987 As- 
chwandcn 201 ll see the Appendix for a brief introduc¬ 


tion to SOC and the relevant statistical method to be 
used in the present work). The energy dissipation driven 
by magnetic reconnection i s believed to be in an SOC 
system, such as solar flares (Lu & Hamilton 1991), pos¬ 
sibly the X-ray flares o f y-ray bursts (CTR Bs ) and black 


holes of various scales (Wang & Dai 2013 Wang et al. 


20151. An SOC system experiences scale-free power-law 


distributions of various event parameters, such as the 
total energy (or equivalently the fluence), the peak en¬ 
ergy dissipation rate, the time duration, and the flux (or 
equivalently the CR) of events. These statistical prop¬ 
erties are related to the geometric dimension of the sys¬ 
tem that drives the events, a nd can thus be used to di - 
agnose their physical nature (Aschwanden 2012 2014). 
While the statistical analysis of solar flares suggests a 
spatial dimension S = 3, consistent with the complex 
magn etic structure observed in active regions on the Sun 
(e.g., iPriest &; Forb es 2000), the X -ray flares of GRBs 
favour S = 1 (|Wang ik Dai||2013|). Thus it is argued 
that the X-ray flares occur in the jets of GRBs where 
the magnetic field conf iguration tends to b e poloidal and 
thus o ne-dimensional (|Wang & Dai||2013|). [Wang et al. 


(2015) further claim S = 3 from an analysis of the X-ray 
flares of Sgr A*, Swift J1644+57, and M 87. However, 
this analysis is very crude in the sense that the correla¬ 
tion among different bins of the cumulative distribution 
are not properly addressed, that artificial cuts of the fit¬ 
ting range of the data are adopted, and that the used 


flare sample is stron 


gly biased and incomplete at the de- 


2004). All these issues are handled 


tection limit (Wang 
in the present work. 

The rest of this paper is organized as follows. The 
observational data and the statistical method adopted 
are described in §2. We present the fitting results in §3, 
and discuss the physical implication of the results in §4. 
Finally we conclude this work in §5. 


2. DATA AND METHODOLOGY 
2.1. Observations 

We use the data of the 2012 Chandra XVP campaign, 
which consists of 38 Chandra ACIS-S/HETGS observa¬ 
tions of Sgr A* with the total exposure time of 3 million 
second (Ms) between February 6 and October 29 in 2012. 


The 2—8 keV light curve is extracted in 300 s bins includ¬ 
ing both photons of the 0 th order and the ±l st diffraction 
orders. The 0 th order events are extracted from a circle 
region with a radius of 1".25, while the ±l st orde r even ts 


are extracted from a 2".5 rectangular region (Neilsen et 


al. 2013). For more details of the XVP campaign and 


the data reduction, we refer the readers to |Neilsen et al.| 

pl3i. 


2.2. Synthetic Light Curve 

The X-ray light curve can be decomposed into the qui¬ 
escent and flare components. The quiescent emission is 
assumed to be steady with a count rate (CR) r. We 
model the flare component as a sum of Gaussian func- 
tion;0 Then the model light curve is 




Ei 


—( a/ 27TCTj 


exp 


1 / 1 - Hi 

2 V a i 


2 


+ r, 


(i) 


where Ei, at, and /ij are the fluence, width and peak 
location of the fth flare, [k] is the integer part of tt and 
it presents the total number of flares. We use a Monte 
Carlo method to generate the model light curve, which 
accounts for the Possion fluctuations of photon counting. 
The differential flare fluence function is assumed to be a 
power-law form, 


, F,~ a B 

dN/dE oc ^ q 


otherwise 


( 2 ) 


where the lower limit of the fluence E m ; n = 1 cts is 
well b elow the detectable flare fluence in INeilsen et al.l 
(2013), and the upper limit .E max = 1000 cts is slightly 
larger than that of the brightest flare (Nowak et al. 2012). 
Given «e and n, the fluences of the model flares can be 
easily sampled. 


As shown in Neilsen et al. (2013), there is a clear cor¬ 
relation between the flare fluence E and duration T (de¬ 
fined as 4 times of the Gaussian width a). We char¬ 
acterize this correlation with a log-normal distribution 
N(log ai, <jet), where ctet is the Gaussian dispersion and 
a power-law relation between <j,; and Ej as 


a i = A [ —— 
1 100 


(3) 


where A is the normalization constant and «et is the 
power-law indes^J The intrinsic scattering of the flare 
widths around the power-law relation is estimated to be 
(Jet = 0.25 for the detected sample (Yuan et al. 2015, 
in preparation). The peak time m of each flare is gener¬ 
ated randomly in the observational periods. The discrete 
light curve is then generated in 300-s bins, resulting in 
M = 9964 bins for the total exposure of the data. The 
average count rates are calculated in each bin of the sim¬ 
ulated light curve. In total we will have five parameters 


2 This modelling is used as the l st -order approximation to ac¬ 
count for the correlation among photons from individual flares, 
although some bright ones do show signific ant asymmetric shapes 
i Baganoff et al. 200l| [Nowak et al.||2012) , and/or substructures 

(hSarriere et al. pmf -— 

15 In this work we use ax to represent the power-law index of 
the distribution of variable X (specifically dN/dX oc X ~ ax ), and 
qxy to represent the correlation power-law index between X and 
Y (specifically Y oc X QXY ). 
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to describe the synthetic light curve, namely the Poisson 
background r, the power-law index of the fluence distri¬ 
bution aEi the total flare number k, the duration nor¬ 
malization A , and the fluence-duration correlation slope 


Qet- 


Finally we correct the pilcup effect of the simulated 
light curve by assuming the l st -order count rates a con- 


20l2] see also 

Neilsen et al. .2013 

20151. For all the sim- 

ulated CRs, tJ 

le pileup correction is less than 20%. 


Fig- 0 shows a representative simulated light curve 
(blue), compared with the observational one (red) with 
the observing gaps removed. 


2.3. Statistical Comparison 

The synthetic light curve cannot be directly compared 
to the observational one due to the lack of bin-to-bin 
correspondence between the two. We need a statistical 
way to make the comparison. A “first” order statistics 
is the CR distribution, which reveals the overall magni¬ 
tude of the va riabilities. However , as mentioned above 
(also noted in Neilsen et al. 2015), the CR distribution 
contains no information about the flare flux correlation 
among adjacent bins. This correlation can be indepen¬ 
dently accounted for by the use of the auto-correlation 


-— ? -_- 

Simonetti et al. 

1985 

Emmanoulopoulos et al. 

2010) 


mine the model parameters. 


2.3.1. Count Rate Distribution 

With the 300-s bins both for the simulated and ob¬ 
served data, the binning results in M = 9964 data 
points for the total exposure window. We construct the 
CR distribution of th e data l ogar ithmically and follow 
Knuth (2006; see also Witzel et al. 2012) to determine 
the optimal bin number of the histogram. Considering 
the histogram as a piecewise-constant model, the rela¬ 
tive logarithmic posterior probability (RLPP) for the bin 
number m, given M data points, is 


( 7TI \ 

— J — m log] 

m / -j \ 

-lo g r(M+-) + 5>gr(ni + *J , (4) 


where ni is the histogram value of the Zth bin, and F(:r) 
is the Gamma function. Then the optimal number of 
bins m can be derived through maximizing the above 
RLPP. For the observational light curve of Sgr A*, we 
find m ~ 14. 

With the above binning, we can now construct y 2 
statistics for the CR distribution analysis: 

XCR = (**obs Wim) C (n 0 bs Wim); (6) 

where the symbol T stands for the transpose of the ma¬ 
trix (vector), n 0 b s and n s i m are the vectors of the CR 
histogram values of the observed and simulated data, re¬ 
spectively. The covariance matrix C is calculated accord¬ 
ing to N Monte Carlo realizations 

i N 

C(b j) = ~ ni ) ( n i “ %) > ( 6 ) 

1 fc= i 


where n* is the mean value of the *th bin of N real¬ 
izations. The covariance matrix is involved here to ac¬ 
count for the error correlat ions among different bins due 
to the correlated flare CRs ( Brockwell fe Davis|2009[ and 
see d iscussions of the error estimation in |lNorberg et ah] 
2009). If there is no correlation between two histogram 


bins, C(i, j) = Sijaf, Eq. © is reduced to the standard 
definition of y 2 - The exact value of IV, as long as it is 
sufficiently large, has little effect on the calculation of the 
covariance matrix. In this work we adopt N = 500. 


2.3.2. Structure Function 


Following Emmanoulopoulos et al. (2010) we define the 
normalized SF as 


V(t) = 


[R(t + r) - R(t)] A 


S 2 


(7) 


where R(t) is the CR as a function of time t, t is the 
time lag, () represents the average over the whole time 
range of the light curve, and S 2 is the variance of the CR 


S ' 2 = 


R(t) - R(t) 


( 8 ) 


This normalized SF is equivalent to the auto-correlation 
function as long as the time series is stati onary, i.e^ 


R(t) = R{t + r) = const. (Emmanoulopoulos et al. 


2010). The SF is especially suitable for the analysis of 


the data that are unequally sampled with large observa¬ 
tional gaps, as is the case here. 

Similar to the CR distribution, we build a y 2 statis¬ 
tics to compare the observed and simulated SFs quanti¬ 
tatively. Apart from the correlation among various SFs 
which needs to be taken into account more seriously com¬ 
pared with the CR distribution, another issue for the 
standard definition of the y 2 is that the random vari- 
a ble (the SF here) needs to be Gaussian. As shown in 


Emmanoulopoulos et al. (2010), the logarithm of the SF 
instead of the SF itself is approximately Gaussiaifi] The 
y 2 statistics for the SF is thus defined as 


Xsf = (log V obs - logV sim ) T C SF 1 (logV obs - log V sim ), 

(9) 

where the covariance matrix Csf can be calculated the 
same way as Eq. ©.• The number of the SF bins for 
the fitting procedure is chosen to be the same as that of 
the CR bins, i.e., m = 14. Note that the SF becomes a 
featureless plateau when the time lag is larger than the 
largest time scale of the flares. Correspondingly, the SF 
for the 14 bins are calculated between 300 s and 3x 10 4 s|3 


4 Although |Emmanoulopoulos et al. | (|2010| show that the loga¬ 
rithm of the Sr is only approximately Gaussian below a false break 
time scale related with its power spectrum, our Monte Carlo simu¬ 
lations suggest that log V (r) is indeed distributed roughly normally 
for the time scales considered here in spite of its ignorance of the 
location of the false break. 

5 Firstly, the SFs for the larger time delay r will show remarkable 
fluctuations, and thus can be impossibly used for the fitting. Sec¬ 
ondly, we calculate the SFs with the observing gaps removed, which 
could introduce an artificial dip in the SFs when r > 8x 10 4 s, 
corresponding to the typical time scale of the observing windows 
for the XVP data. However, the observing gaps show no significant 
effect on the SFs with r<3x 10 4 s, the fitting range considered 
here. 
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Fig. 1 . — Upper panel: The light curve (with the observing gaps removed) from the XVP observations (red), compared with 
a representative one from the simulations (blue) using the best-fit parameters of model I shown in Table [I] (see below). Lower 
panel: Zoomed-in light curve of a time window indicated by the short horizontal bar in the upper panel. 


By minimizing y * 2 = + x| F , we test our model and 

constrain the model parameters. 


2.3.3. Model Fitting 

We minimize \ 2 = Xcr + X'sf> using the MCMC 
method. It maps out the full posterior probability distri¬ 
butions and correlations of the model parameters. The 
MCMC code is adapte d from the public CosmoMC code 


rithm ([Metropolis et al. 1953 


(Lewis fe Bridle 2002). The Metropolis-Ha stings algo- 


and-accept process m which 


Hastings 1970), a propose- 
;he acceptance or rejection 


of a proposed point in the parameter space depends on 

the probability ratio between this point and the previous 
one, is adopted to generate the Markov chains. 


3. RESULTS 
3.1. Joint Fitting 

The parameter ctet hr the fittings is chosen to be 0.25 
according to Yuan et al. (2015). Fig. [2] shows the best-fit 
results and la bands of the CR (upper panel) and the SF 
(lower panel), compared with the observations. The CR 
distribution exhibits an exponential distribution, which 
reflects the Possion background, and a power-law tail 
which represents the flares. The SF increases with r, and 
reaches a plateau of r ~ 5 x 10 3 s, which corresponds to 
the largest duration of the flares. The profile with a ris¬ 
ing trend clearly shows the time scale of the flares of the 
light curve. When the time lag r is much longer than 
the largest time scale of the variabilities, the normalized 


SF asymptotically approaches 2, i.e., the plateau in Fig. 
EJ It should be noted that the SF can be biased when 
applied to irregularly sampled light curves. In particu¬ 
lar, the SF can show spurious breaks at At ~ 10 _1 AT, 
w here A T is the tot al exposu re time of the light curve 
(Emmanoulopoulos et al. 2010). However, the break seen 
here is at a much shorter time scale (~ 10~ 3 AT), and 
consisten t with break time scales seen previous ly in both 
the NIR ( MeyerjTjdj[200 9^ Witzel et al. ||2Q L2] ) and sub- 


millimeter ( |h)exter et al. 


that the turnover time is 


21)14 b Both of them suggest 
ikeiy an intrinsic characteristic 


time scale for the X-ray flares. 

Fig. [3] shows the one and two dimensional (ID and 
2D) probability distributions of the parameters. The 
mean values and the ler uncertainties of the parame¬ 
ters are listed in Table |T| The fitting results suggest 
a fluence distribution power-law index ag = 1.65 ± 0.17, 
and the less constrained fluence-duration correlation in¬ 
dex q:et < 0.55 (95% confidence limit). The correla¬ 
tions between some of the parameters can be clearly seen 
in Fig. [3} There is a strong correlation between the total 
number of flares k and the fluence distribution index qe- 
It is easy to understand that a harder fluence distribution 
will naturally correspond to a smaller number of flares in 
order not to over-produce the number of photons. Also 
the anti-correlation between (or k) and r is again due 
to the constraint of the total number of photon counts. 

To judge the “goodness of fit” of the best-fit model, 
we use the bootstrapping method to estimate the confi- 
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Fig. 2. — Upper panel: The CR distribution of the X-ray 
light curves of Sgr A*. The red circles are for the observed 
light curve, and blue dots are the best-fit model results. The 
bands show the fitting la ranges of the CR distribution for 
the best-fit model according to the MCMC chains. Lower 
panel: The SF for the observed data and the simulations for 
the same parameter set. 


TABLE 1 

Best-fit model parameters with the MCMC method 


r (10 3 cts s - 1 ) 

QE log(K) 

log(A) 

OET xi 

5.90 ±0.14 

1.65 ±0.17 2.41 ±0.22 

2.99 ±0.12 

< 0.55 0.9 


Note. — The errors are all lcr (or 68%) confidence limit 
except for oet (for the one tail 95% confidence). 


dence level. Based on the best-fit model parameter set 
as shown in Table [ll we generate 2000 realizations of the 
light curves and calculate the y 2 values for each real¬ 
ization using the same method described in §2.3. The 
distribution of the \ 2 values is shown in Fig. [4] The 
number fraction of realizations with x 2 values smaller 
than that of the observational one, x 2 bs = 21, is esti¬ 
mated to be 89.7%, indicating that the observed light 
curve is well characterized by the best-fit model. 

Our preliminary re-analysis of the light curve indicates 
that the existing sample of flares may be biased against 
the detection of long duration ones. This bias can sig¬ 
nificantly affect the measurement of <tet and potentially 
other parameters as well. We have thus tested this lat¬ 


ter sensitivity by fixing <tet to different values and find 
that the estimates of other parameters are not signifi¬ 
cantly altered (well within their uncertainties). The anal¬ 
ysis, however, may be sensitive to the assumed specific 
shape of individual flares. We will address these poten¬ 
tial higher-order complications in a follow-up paper. 


3.2. Comparison with previous works 


Neilsen et al. (2013) studied the statistical properties 
of the iiuences, peak rates, durations, and luminosities, 
of the 39 flares detected during the 2012 Chandra XVP. 
It is found that the distributions of the durations and lu¬ 
minosities can be well described by power-law functions, 
with indices q:t = 0.9 ± 0.2 and o:l = 1.9(tg'3, respec¬ 
tively. The power-law fittings to the fluence and peak 
rate distributions give cce = 1.5 ± 0.2 and op = 1-9 +q( 4. 
However, the fittings can be significantly im proved by as¬ 
suming cutoff power-law functions (Neilsen et al.] |2013 ). 
Our result of qr are consistent with that derived in 
Neilsen et al. (2013|, even though only the detected sam- 
ple was employed in the latter. Based on the Gaussian 
flare profile assumption, the slopes of the duration and 
peak rate distributions can be derived according to our 
fitting results®] which are ap > 2.1 and 1.7 < ap < 2.4, 
respectively. TRe p eak ra te dis tribution is consistent 
with that of Neilsen et al. (2013), but the duration dis¬ 
tribution is different. As mentioned above, the fluence- 
duration correlation derived in our work is also differ¬ 
ent from that of the detected sample. Possible reasons 
for these discrepancies include the i ncompleteness of the 
sample, and the profiles of the flares. Neilsen et al. (20151 
showed that the flux distribution can be decomposed 
into two components, a steady Poisson background and a 
variable flare component. The X-ray flares excess follows 
a power-law process with of = 1.92to)o2 an d the Poisson 
background rate is (5.24±0.08) x 10 -3 cts s -1 . The back¬ 
ground rate derived here (5.90 ± 0.14) x 10 -3 cts s _1 is 
slightly larger than that derived in Neilsen et al. (2015), 
possibly due to the different methods adopted. Although 
not shown explicitly, our result of Fig. [2] gives roughly 
qf ~ 2. 2 which is consisten t with that of |Neilsen et al.| 
(|2015|). | Wang et al.| (120151) a lso analyzed the detected 
sample of Neilsen et al. ~[ 2013), with different selection of 
the fitting data ranges, andgave that a E = 1-8 ±0.6 and 
q:t = 1-9 ± O.fQ Within the large uncertainties, their 
results are consistent with what we obtain. 


4. IMPLICATIONS ON THE NATURE OF THE 
FLARES 

4.1. Confronting the Statistical Results with SOC 
theory 

We find that both the flare fluence distribution and 
the correlation between the duration and the fluence can 
be described by power-laws, which is consistent with the 

6 For the duration distribution, since T oc E a ET , we have = 

al 

H • oc E~ a E • E-“et+i oc T -i-(c*E-l)/a E T i an d thus o T = 
1 + (q;e —1)/o:et- For the peak count rate distribution, P oc E/T oc 
E/E a ET = E 1 —^ ET , then ap = (cke — Q: et)/(1 — »et)- 

7 Concerning their analysis of Sgr A*, however, the flares sample 
isn’t corrected for the incompleteness, as well as the correlation of 
the uncertainties among different bins of the cumulative distribu¬ 
tions are not properly addressed in the fittings. 
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Fig. 3.— ID and 2D probability distributions of the fitting parameters. The contours in the 2D plots are for confidence levels 
68% and 95% from inside to outside, respectively. 



Fig. 4.— The distribution of y 2 values of 2000 realizations 
of the light curves with the best-fit model parameters. The 
shaded region corresponds to a percentile of 89.7% for y 2 < 

Xobs ■ 


prediction of the frac tal-diffusive SOC theory (Aschwan- 


denpOlT 2012 2014). The SOC theory further predicts 


that the power-law indices will depend on the Euclidean 
space dimension S of the system to produce flares. The 
predicted power-law indices are a ~ 1.5 and ss 0.5, 
respectively, for S = 3, the classical diffusion parameter 
fl = 1 and the mean fractal dimension D$ = (1 + S)/2 
(see the Appendix for more details). These values may 
vary in a range for different assumed values of ft and D$. 
For example, «e can range from 1.4 to 1.7 for the fractal 
dimension 1 < Ds =3 < 3. Our results of «e and «et are 
actually in good agreement with the theoretical expecta¬ 


tion with S = 3. The predicted indices of the duration 
and peak rate distributions are a^! 1 ss 2.0 and ap ~ 1.7. 
As a comparison, our induced values are ap > 2.1 and 
1.7 < ap < 2.4, respectively. There are potential dis¬ 
crepancies of these two distributions. Since the deter¬ 
mination of the peak rate depends on the assumption of 
the flare profile, as well as the precise measurement of 
the flare duration, there should be relatively large uncer¬ 
tain ty of the peak CR distribution. As shown in jNowak| 


et al. (20121, the profile, at least for the bright ones, is 
indeed asymmetric rather than Gaussian. Therefore, the 
integral property (fluence) should be more reliably mea¬ 
sured and more suitable to be used to compare with the 
theoretical model expectation. 


4.2. Episodic ejection of plasma blobs as origin of 
flares ? 

Two main conclusions can thus be obtained from the 
above analysis. First, the power-law distributions of the 
fluences, durations, and their correlation, suggest that 
the flares of Sgr A* can be explained in the fractal- 
diffusive SOC framework. Second, the inferred space di¬ 
mension responsible for the flares is S = 3. Both results 
are similar to the solar flares, which thus implies that 
the X-ray flares of Sgr A* are likely driven by the similar 
physical mechanism as that of the solar flares. By anal¬ 
ogy wi th the coro nal mass ejections and their solar flares, 
|Yuan et al.| (2009) have proposed a magnetic reconnec¬ 
tion model for the episodic ejections and the production 
of flaresj^] from the accretion flow of Sgr A*. In the fol¬ 
lowing, we briefly review the key points of the model. 

The structure of the accretion flow of Sgr A* is quite 


8 The flares and episodic ejections are physically associated with 
each other, both for the Sun and the black hole accretion flow. 
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similar to the atmosphere of the Sun, i.e., a dense 
disk enveloped by a tenuous corona, as shown by the 


m agneto- hydrodynamic simulations (e.g., Fig. 4 in Yuan 


& Narayan 20141. The magnetic loops emerge into the 


disk corona due to the Parker instability. The configu¬ 
ration of the coronal magnetic held emerging from the 
accretion flow is controlled by the convective turbulence 
motion of the plasma in the disk. Since the foot points of 
the Held lines are anchored in the accretion flow which is 
turbulent and differentially rotating, a swamp of small- 
scale magnetic reconnection sets in, which redistributes 
the helicity and stores most of it in a flux rope. The tur¬ 
bulent processes in the accretion flow thus continuously 
build up magnetic stress and helicity. When a thresh¬ 
old is reached, e.g., when the current density inside the 
current sheet below the flux rope is strong enough, micro¬ 
scopic instabilities such as the two-streaming instability 
would be triggered, resulting in anomalous resistivity and 


f ast m agnetic reconnection inside the current sheet (Chen 


2011). The equilibrium of the flux rope breaks down 


with accompanied dissipation of magnetic energy in a 
catastrophic manner, powering the observed Hare. This 
is the so-called much-awaited SOC state which means 
that a small perturbation, owing to turbulence motion 
in accretion flow and subsequent magnetic reconnection, 
will trigger an avalanche-like chain reaction of any size 
in the system once the system self-organizes to a critical 
state. 


5. CONCLUSIONS AND DISCUSSION 

We present a statistical analysis to the Chandra 2 — 8 
keV X-ray light curve of Sgr A* from the 3 Ms XVP cam¬ 
paign. A Monte Carlo simulation method is adopted to 
generate the model light curves taking into account for 
the Possion fluctuations of both the photon counts and 
the number of flares. Then we fit the CR distribution 
and the SF of the light curves jointly to constrain the 
model parameters, via an MCMC method. We find that 
the X-ray emission of Sgr A* can be well modelled by two 
distinct components: a steady quiescent emission around 
the Bondi radius with a CR of (5.90±0.14) x 10~ 3 cts s -1 , 
and a flaring emission with a power-law fluence distri¬ 
bution dN/dE oc E~ aE with cte = 1-65 ± 0.17. The 
duration-fluence correlation can also be modelled by a 
power-law form T oc E“ ET with «et < 0.55. 

These statistical properties are consistent with the the¬ 
oretical predications of the SOC system with the Euclid 
spatial di mens io n S = 3, same as that for solar flares (As- 
chwanden 2012, 2014). Our analysis, therefore, indicates 
that the X-ray flares of Sgr A* are possibly driven by 
the same physical mechanism as that of t he s olar flares, 
i.e., magnetic reconnection (Shibata & Magara 2011). 
This idea is further supported by the recent development 
of magneto-hydrodynamic simulations which lead to the 
consensus that the accretion flow in Sgr A* is e nveloped 
by a tenuous c orona above the dense disc (e.g., Yuan & 
Narayan 20141, similar to the atmosphere of the Sun. 
The three dimensional geometry of the energy dissipa¬ 
tion domain further suggests that the X-ray flaring of 
Sgr A* occurs in the surface of the accretion flow due to 
the less-ordered magnetic field structure embedded in it 
compared with tha t in the relativis tic astrophysical jets 
(e.g., the GRBs; Wang & Dai| |2013 ). 

The flares of Sgr A have also been detected in the 


NIR band (e.g., Genzel et al.[2003 GhezetaL 2004). 


Dodds-Eden et al? 1 20111 and|Witzel et al. 2012) have 
analysed the MR flare data and found that the NIR 
flux distribution can also be described by a power-law 
form. However , the power-law indices obtained, which 
are —2.1 ± 0.6 (Dodds-Eden et al. 12011 ) and —4.2 ± 0.1 
(Witzel et al. 2012), are different from each other. The 
cause for this discrepancy between these two results is 
unclear, and may be related to the subtractions of the 
stellar light. The forme r is roughly consist ent with that 
of the X-ray flares (e.g.,|Neilsen et al.|2015| as well as our 
results in the context of the SOC framework with S = 3). 
If this is the case, then the NIR flares are expected to 
have the same physical origin as the X-ray flares, which 
is also supported by the results that NIR and X-ray flares 
occur simultaneously when the re are coordinated obser - 
vations at two wavebands (e.g., Dodds-Eden et al.|2009 ). 
Otherwise, it may suggest different radiation mechanisms 
between the NIR and X-ray flares which lead_to the dif¬ 
ferences in the fractal dimension ( Aschwanden 2011) or 
even different physical origins (Chan et ai.|2015|T 
By analogy with solar flares, the multi-wave band flares 
of Sgr A* in X-ray, probably NIR, and less prominent 
sub-millimeter and radio, are likely associated with the 
ejection of plasmoids, both of which are the radiative 


Yusef-Zadeh et al. 

2006 

2008) 

Brinkerink et al. [2015 


gins of the flares and episodic jets in various black holes 
in general. 


Mineshige et al. ( 1994b; see also Negoro et al. 1995 


Takeuchi et al. 1995) studied the statistical properties 


of the X-ray fluctuations in Cygnus X-l, a black hole 
X-ray binary. They found that the X-ray fluctuations 
showed much steeper (or even exponential) distributions 
compared with those of solar flares. They argued that 
the fluctuations should also follow an SOC process. A 
mass diffusion process was invoked in their works to ex¬ 
plain the discrepancy of the distributions. But the X-ray 
fluctuations they discussed likely have a different physi¬ 
cal origin compared with what we have discussed in the 
present paper. Thus the trigger mechanism of SOC is 
also different. 
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APPENDIX 


SELF-ORGANIZED CRITICALITY SYSTEM 


The theoretical concept of SOC, proposed first by 
class of nonlinear dissipative dynamical systems whic 


Katz (1986) and Bak et al. (1987) independently, describes a 
1 are slowly and continuously driven toward a critical value of 


an instability thres hold, producing scale-free, fractal-diffusive, and intermitten t avalanches fCh arbonneau et al.||2001 
Aschwanden 2011). The classical example of the SOC system is a sandpile (|Bak et al.||l987h . If one continuously 


2011|). The classical example of the SOC system is a sandpile (Bak et al. 1987k If one continuously 
drops the sand grains (driving force) to the same place, a conical sandpile will grow m a steady way until the surface 
shape reaches a critical slope, beyond which further addition of sand rapidly leads to catastrophic avalanches of the 
sandpile. The critical slope is primarily determined by the friction threshold between adjacent grains. Such a threshold 
is crucial since it allows the existence of multiple metastable states across the system. The continuous slow addition 
of sand will produce small or large avalanches with sizes independent of the input rate of sand. The critical behavior 
of the sandpile appears naturally as a consequence of the slow addtio n of sand without fine tuning of the addition 


rate. It is in this sense that the system is self-organized to criticality (Charbonneau et al. 2001 Aschwanden 


2011 ). 


The continuous energy input (addition of sand) and nonlinear energy dissipation (avalanches due to the complicated 
interactions between the colliding sand grains) are two key points to determine whether a system is in an SOC state. 

A universal property of the SOC system is that there is no preferred scale of the release of energy. The scale- 
free power-law distributions of various event parameters, which actually become the hallmarks of SOC systems, are 
expected. The ubiquitousness of power-law distributions in a wide range of physical systems, e.g., earthquakes, cloud 
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formation, solar flares, and widely in a strophysics, sug gests that SOC is a common la w of the nat ure. The SOC theory 
has b een widely applied in geophysic s (Turcotto 1999), solar physics (e.g.,[C harbonnea u et a l.||2001|), and astrophysics 
(e.g., Mineshige et al.||1994a Kawaguchi et al.||2UUU ' Aschwanden||2011 Kunjaya et al. 2011| ). 

The slope of the power-law form depends on how the subsystem sel f-organizes to the c ritical state. An analytic 
description was provided based on a fractal-diffusive avalanche model (Aschwanden 2012). The SOC state ensures 
that the entire system is close to the instability threshold, and the avalanches can develop in any direction once 
trig gered. The propagation of the unstable nodes can thus be modelled with a random walk in an 5-dimensional space 
(Aschwanden 2012J). A statistical correlation between the spatial length scale L and time duration of the avalanche 
T for a diffusive random walk is predicted to be L oc T^ 3 / 2 , where /3 is a diffusion parameter [j3 = 1 corresponds to 
the classical diffusion). With the scale-free probability conjection dN/dL oc L _s , which means that the critical states 
are homegeneously distributed across the entire system, the occurrence frequency distri butions of flare duration, flux, 
peak flux, and total energy can be obtained in the fractal-diffusive avalanche model (Aschwanden||2012). The indices 
of the total energy (fluence) off, du rati o n off , peak luminosity (or peak CR in this work) ap 1 , and flux (or CR) ap 
distributions are derived as (Aschwanden|j2012 2014) 


th_i 

a E —1 

_ th _ i 

Oi ’p — 1 

_ th_i 

Qlp — _L 

_ th_i 

CUp — 1 


5-1 

+ A3 + 2//3’ 
, {S-l)p 
2 

5 — 1 

+ ^T’ 

5-1 

H-^—5 


(Al) 

(A2) 

(A3) 

(A4) 


where 5 is the Euclidean space dimension, Dg is the fractal Hausdorff dimension which lies between 1 and 5. The 
fractal-diffusive SOC theory also predicts the s caling relations between v arious flare parameters. The correlation slope 
between the total energy E and duration T is (Aschwanden||2012 2014) 


th _ 

a ET — 


DsP + 2 


(A5) 


For the classical diffusion (/3 = 1) and an estimated mean fractal dimension of Dg ss (5 + l)/2, we have ap = 1.5, 
c4 h = 2.0, a)! 1 = 1.7, ap = 2.0, and a^ h T = 0.5 for 5 = 3. 



































